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ABSTRACT 

We re-examine the effects of redsliift space distortion in all-sky galaxy redshift surveys 
in the formalism of spherical harmonics. This natural decomposition of the density field 
into radial and angular eigenfunctions of the Laplacian operator complements both the 
spherical symmetry of the survey geometry and the dynamical basis for the redshift 
distortion. 

Within this framework we show how one can treat both the large-scale linear 
effects, and the small-scale nonlinear clustering, exactly to first order. We show, contrary 
to earlier claims, that the redshifted density field is no longer homogeneous, as well as 
being anisotropic. 

We construct a likelihood function for each mode in the decomposition, based 
on the random Gaussian field hypothesis. This function finds its maximum when the 
underlying field is homogeneous and isotropic, and has the correct amplitude at each 
mode. As the level of distortion in the field is governed by the cosmological density 
" ' parameter, via j3 = Slg'®/&, where b is the galaxy bias parameter, strong limits can be 

placed on cosmological models, 
r The method also allows in principle a determination of the power spectrum of per- 

O^l turbations, requiring no assumptions beyond that of linear theory. The method therefore 

offers significant advantages over Fourier techniques when dealing with all-sky surveys. 

We apply our likelihood analysis to both simulated data, and real data, using the 
IRAS 1.2-Jy galaxy catalogue, for which we find a maximum likelihood /3 ~ 1.1 ± 0.3, 
and a real-space fiuctuation amplitude corresponding to (T8_iras = 0.68 ± 0.05. The l-a 
errors should be treated cautiously and are discussed in the paper. 
i-C We also relax the gravitational instability assumption, to find a more general de- 

termination of the velocity power spectrum required to reconcile the anisotropic redshift 
^ space map with the assumed isotropic real-space map. 
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1 INTRODUCTION 

An accurate three-dimensional map of the galaxy distribu- 
tion would be enormously valuable for cosmology, but the 
lack of an accurate estimator of distance precludes this. De- 
spite best efforts, the errors in the TuUy-Fisher and Faber- 
Jackson distance indicators mean that the best one can do 
is to use the recession velocity v, or redshift, of a galaxy, and 
to assign its distance s by using the Hubble expansion law 
V = Hos, where Ho is the Hubble parameter. The resulting 
'redshift space' map is of course not perfect, because galaxies 
are not necessarily moving precisely with the idealised ex- 
pansion of a homogeneous universe. Density inhomogeneities 
have associated with them peculiar velocities, which intro- 
duce a radial distortion between the true (real-space) map 
and the redshift space map. This distortion is, of course, in- 
convenient: for example, extracting the power spectrum of 



density fluctuations from a galaxy redshift survey becomes 
a non-trivial task. However, the distortion itself can be ex- 
ploited, as the magnitude of the distortion is dependent on 
the density parameter of the Universe Qo, if fluctuations 
grow by gravitational instability. Constraining Qo is one of 
the main aims of this paper. 

In the most general terms, we can divide the meth- 
ods for constraining Qo with peculiar velocities into two 
catagories which we shall label dynamical and statistical. 
In the former catagory perturbations in the density field are 
related to the peculiar velocity field through the continuity 
equation. Notable examples of this type are dipole studies 
(e.g. Lynden-Bell et al. 1989, Rowan- Robinson et al. 1990, 
Strauss et al. 1992), comparing the absolute motion of the 
local group to the gravitational field implied by redshift- 
surveys; the POTENT analysis (Bertschinger & Dekel 1989, 
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Dekel et al. 1993), making a point-by-point comparison of 
the reconstructed density field with redshift-maps; and the 
inverse problem of comparing the predicted to measured ra- 
dial, relative peculiar velocities from redshift space maps 
(Kaiser et al. 1991). The analysis is complicated by the need 
to correct the density field for redshift-distortions, and in the 
latter two cases, the need to know what the radial peculiar 
velocities actually are. Given the uncertainties in redshift- 
independent distance indicators, and an ongoing debate over 
the effects of Malmquist bias, it is valuable to have a method 
which is independent of peculiar velocity measurements. 

The need to measure peculiar velocities can be over- 
come within a statistical framework, and hence offers a 
powerful alternative to the dynamical methodology. The 
redshift-distortion introduces spurious terms into the sta- 
tistical properties of the galaxy density field, which can be 
measured and used to constrain Qq. In this paper we shall 
develop the statistical approach. 

Coherent galaxy motion has the effect of locally trans- 
forming the real space coordinate system: 

(r, 0, <f>) (s, e, <t>) = (r, 0, <f>) + (u, 0, 0), (1) 

where radial distances are estimated from Doppler redshifts. 
Here u = v(r).r/_ffo is the radially projected galaxy peculiar 
displacement. An order of magnitude estimate suggests that 
the fractional change in the linear density field due to this 
transformation will be ~ /? = Qg^/fc, the ratio of linear 
growth to galaxy bias parameter . Hence, by measuring this 
change in the density field we can hope to constrain /?. 

In a seminal paper. Kaiser (1987; see also Lilje & Efs- 
tathiou 1989, McGill 1990) has given an approximate anal- 
ysis of the distortion effect, which has been the basis for 
subsequent work. The Kaiser analysis consists of decompos- 
ing the density field into plane-waves, arriving at the linear 
expression 

S'(kln) = (l+l3lI^)S(klr,), (2) 

in agreement with our approximate result, where fi is the 
cosine of the angle between the plane-wave and the line 
of sight. For high Qo, as expected from inflation, and 
low galaxy bias, the projection into redshift space is non- 
negligible - the redshift space density differs from the the 
real-space density field by factors of order unity. 

The distortion equation is derived assuming the gravita- 
tional instability hypothesis - now on a firm footing since the 
discovery of small initial fluctuations implied by the COBE- 
DMR observations, but can be generalised to models with 
non-gravitational velocities, provided that they are curl-free. 
However, the Kaiser equation can really only be applied uni- 
versally if the survey subtends a small angle, so that the 

* Note that all methods which relate deviations from uniform 
density to deviations from uniform expansion are bound (at least 
within linear theory) to give answers for Qq which are uncertain 
because of ignorance of how the measured (number) overdensity 
of galaxies is related to the mass overdensity Sp. The fractional 
number overdensity Sn/n is assumed in simple analyses to be 
bSp/p, where 6 is a uniform bias parameter. In our analysis b is 
defined in terms of the power spectrum. This is a weaker condi- 
tion than relating number densities, as only the long wavelength 
modes used in our analysis are assumed to have a simple ampli- 
tude enhancement. 



radial distortion can be approximated by a distortion along 
one axis. For a large-angle survey, one may be able implicitly 
to split the survey into small, independent volumes, but this 
will surely not be possible if the wavelength of interest sub- 
tends a large angle. Herein lies the problem of the Fourier 
approach: one wants to analyse short wavelengths, so the 
small-angle approximation may hold, but it is these short 
wavelengths which are difficult to analyse because they are 
nonlinear. With current survey depths, it is by no means 
obvious if there are any wavenumbers at all which are both 
linear and small-angle (cf Cole, Fisher & Weinberg 1994). 

Further problems arise due to the choice of decompo- 
sition into Cartesian Fourier modes, a common choice for 
calculating the power spectrum from the galaxy distribution 
(we shall have more to say about statistics below). As the 
distortion is purely radial, we expect mixing of the Fourier 
modes (Zaroubi & Hoffman 1994). This will also arise as 
a result of finite and nonuniform sampling of the galaxy 
density field. Generally these effects can be reduced to in- 
dependent radial and angular selection, the latter arising 
from incomplete sky-coverage. The redshift space distortion 
only interferes with the radial component, producing den- 
sity distortions and surface effects. Finally, we would expect 
a purely radial distortion to destroy the statistical homo- 
geneity and isotropy of the density field, whereas the Kaiser 
equation, while destroying isotropy, preserves homogeneity. 

In this paper we re-examine the effect of redshift- 
distortion on a finite, inhomogeneously-sampled galaxy dis- 
tribution using a decomposition into spherical Bessel func- 
tions and spherical harmonics, a formalism which was sug- 
gested for a somewhat different problem by Binney & Quinn 
(1991), and whose use for analysing galaxy catalogues was 
noted by Lahav (1993). The effects of redshift distortions 
in this description have been found independently by Fisher 
et al (1994c) and applied to reconstructing the real-space 
density field. 

The method presented in this paper extends work in 
spherical harmonics in cosmology which goes back to Pee- 
bles (1973), and which includes work on analysis of near 
all-sky catalogues (Fabbri & Natale 1990, Scharf et al. 1992, 
Scharf & Lahav 1993) and velocity fields (Regos & Szalay 
1989). The fact that redshift distortions affect the observed 
harmonics was noted by Scharf & Lahav (1993), and the dis- 
tortion was exploited, without making a full decomposition 
of the field as here, to estimate /? by Fisher et al. (1994a). 

Spherical coordinates are clearly the natural system to 
use, given the nature of the distortion, and a number of 
distinct advantages will become apparent. The two main 
advantages of the approach in this paper are firstly that 
the only approximation which one has explicitly to make is 
that the distortion is linear^ . Since this method is essentially 
a power spectrum approach, it is straightforward to select 
only those modes which have long enough wavelength still 
to be in the linear regime. It therefore avoids the difficulties 
of correlation function methods (cf Hamilton 1992, Davis & 
Peebles 1983, Fisher et al. 1994b) in which the separation 
into linear and nonlinear regimes is not so clear-cut, and 



T A further approximation concerning discreteness of modes is 
implicit in our choice of boundary condition, similar to the usual 
practice of performing Fourier expansions in a cube of finite size. 



Spherical Harmonic Analysis of Redshift Space 3 



which usually involve noisy integrations of the correlation 
function out to large distances. 

The second major advantage is that this formalism of- 
fers a method of estimating /? and the real-space power spec- 
trum itself. This facet is not pursued fully in this paper, but 
it is an important point. One of the major goals of cosmology 
in the next few years will be to measure the power spectrum 
on scales in the large gap between current measurements and 
those probed by microwave background experiments such as 
COBE. With reasonable depth surveys, one will have to con- 
sider the spherical nature of the distortion. 

There are other minor advantages of this formalism: in 
a spherical coordinate system we can separate out the an- 
gular and radial effects and treat them independently, thus 
allowing us to deal with surface effects by imposing suit- 
able boundary conditions. In addition we may calculate the 
change of phase information, as well as amplitude, of the 
harmonic coefficients. And finally, in applications where the 
distribution of objects sample cosmologically significant vol- 
umes, e.g. in the case of measuring structure in quasar sur- 
veys, curvature effects become important. In the framework 
of a spherical harmonic decomposition the radial component 
can be generalised to account for curvature effects, given a 
cosmological model. 

In order to constrain cosmological parameters, we need 
to compare the anisotropic redshift space pattern of galaxy 
clustering with the expected underlying distribution. As we 
have no direct information as to the actual pattern, this 
must be done statistically. To proceed we need an a pri- 
ori assumption for the underlying density distribution. The 
central limit theorem lends weight to the Random Gaus- 
sian Field Hypothesis, either from general considerations of 
the sum of many independently generated perturbations, or 
from inflationary theory. This also has observational sup- 
port from the Rayleigh distributions of powers (Feldman et 
al. f994), topology studies (Gott et al. f986) and studies 
of skewness (Juskiewicz f993). The distribution is then fully 
specified by the (isotropic) power spectrum of perturbations, 
P{k = |k|). This can be either estimated along with /?, from 
the data using a likelihood analysis, or put in by hand to 
find the maximal value of /? for a given model. 

The form of the paper is set out as follows. In Section 3 
and 4 we formulate the description of redshift space density 
perturbations in the language of spherical harmonics. This 
leads to a simple linear matrix relation between the modes of 
the underlying density field and the observed fluctuations in 
galaxy density in redshift-space. The stochastic properties 
of the observed galaxy density distribution are quantified us- 
ing their covariance matrix in Section 4, including the shot 
noise contribution from finite sampling. We also consider 
here the effects of small scale, nonlinear velocity dispersion 
in redshift space. In Section 5 we expand on the stochastic 
properties of the redshift space density field by calculating 
the distribution of density fluctuations. Practical schemes 
for weighting the data are also considered in Section 5. The 
formalism is applied to a likelihood analysis of N-body and 
IRAS 1.2 Jy datasets in Section 6, using a parametric form 
for the underlying power spectrum. We also discuss in Sec- 
tion 6 how other observations, of the microwave background 
radiation, and of the abundances of rich clusters, may be 
combined with our analysis to constrain Qo itself. 



1.1 Outline of the method 

The method of analysis described in the paper consists of 
the following: 

• Expansion of the redshift-space density field in spherical 
harmonics and spherical Bessel functions (equation 13). 

• Conversion of real-space coefficients to redshift space, 
assuming linear theory, and taking into account incomplete 
sky coverage, selection function and real-to-redshift-space 
distortion (equation 18). 

• Assuming an isotropic model power spectrum in real 
space, calculating the covariance matrix of the model coef- 
ficients in redshift space. 

• Using the maximum likelihood technique to estimate 
the distortion parameter /? = OP'^ /h and the underlying 
real-space power spectrum. 



2 PECULIAR VELOCITY FIELDS IN LINEAR 
THEORY 

In a redshift space map, galaxies are placed at a position 
St={s,d,4>), where the distance coordinate s is the recession 
velocity divided by the Hubble constant Ho - In general this is 
not the true distance because the galaxy may have a peculiar 
velocity v. The redshift space position is then related to the 
real-space position r by 



where Vo is the peculiar velocity locally. 

In linear theory, the growing-mode of gravitational per- 
turbations is curl-free (e.g. Peebles 1980), and the peculiar 
velocity may be expressed as the gradient of a velocity po- 
tential: V = — V^*. ^ is related to the mass overdensity 
8p = hpjp by V^* = Hof{Q.o)hp where /(Qo) ~ is the 
logarithmic growth rate with respect to scale factor of the 
growing-mode perturbation. 



3 SPHERICAL HARMONIC FORMALISM 
3.1 Introduction and goals 

In this section we set out the main equations relating the 
density field in redshift-space with the underlying density 
perturbations, and define the notation used throughout this 
paper. The goal is to relate the transform coefficients of a 
weighted redshift space survey to the transform of the actual 
density field in real space. The main result of this analysis is 
that the fields are related by a convolution, cast as a matrix 
relation between real and redshift space density modes. We 
begin by discussing the spherical harmonic formalism. 



3.2 Spherical Harmonics 

Spherical harmonics arise in physics as a result of their 
use as eigenfunctions of the Laplacian operator in systems 
with some spherical symmetry. Hence they can form a com- 
plete set of orthonormal basis functions. Following Binney 
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& Quinn (1991) we expand scalar quantities in spherical 
harmonics via the Bessel-Fourier relation 



A{r) = ClnAlm„Jl{klnr)Ylm{9, (f)), 



A,, 



/ A{v)n{h„r)Yil,{e,4,)d\, 



(4) 



(5) 



where cin is a normalisation factor, chosen to appear sym- 
metrically in the transform and its inverse. ji(y) is the l"* 
order spherical Bessel function, related to ordinary Bessel 
functions by ji{z) = ^Jl^f{2z)J^_^l{z), and Yim{i) is the 
Spherical Harmonic function of order I and m, in the an- 
gular direction r. A number of conventions exist for spher- 
ical harmonics; we adopt that used by Binney & Tremaine 
(1987): 



21 + 1(1- \m\)\ 



47r {l+\m\)\ ' 
(-1 



pl'"l(cose) 



X exp(imi/i) x 



m > 
1 m < 



(6) 



The wavenumbers, kin, are chosen subject to a boundary 
condition on A on the surface of a sphere of large radius Tmax 
(see Section 3.3). Note that both the cin and kin depend on 
I as well as n. 

In the notation of Dirac, we can represent the scalar 
product of the vector basis of configuration space with spher- 
ical harmonic space by (Binney & Quinn 1991) 

(r|lmra) = Clnkln]l(klnT)Ylm(i), (7) 



3.3 Boundary Conditions and Eigenvectors 

As the potential field of the density distribution must satisfy 
Poisson's equation, unique solutions exist only for Dirichlet 
or Neumann conditions at the boundary (Jackson 1975). 
To first order we can avoid surface terms on the boundary 
of the survey by imposing Neumann boundary conditions, 
where the first derivative of the potential field normal to the 
boundary vanishes. Physically this requires that the velocity 
field be zero on the survey boundary, and so there is no 
distortion. 

The surface condition of zero potential gradient implies 
that in spherical harmonic space the wavenumbers kin are 
chosen such that 



d_ 

dr 



Jl(klnT] 



0, 



(8) 



on the boundary. This is a somewhat different assumption 
from that used by Taylor and Rowan-Robinson (1993) when 
reconstructing the density field from redshift space, who set 
the potential to zero on and beyond the boundary. This 
ensured that the survey was treated as an isolated system. 
It also differs from the boundary condition used by Lahav 
(1993) and Fisher et al. (1994c), who required the potential 
to have a continuous derivative at the boundary. In our case, 
the selection function is very small at the boundary, so the 
results should not be sensitive to the choice of boundary 
condition. We have tested that the boundary conditions do 
not affect our results by placing the boundary at 200 and 
250 Mpc, and find consistent results. 

For a finite sphere, we must introduce a normalisation 



constant to preserve normality. This can be calculated us- 
ing the orthogonality properties of the harmonic basis (see 
Appendix C): 



2k' 



3/2 



^^^{[1/4 + fcf„r?^ax - (I + ll2Y]Jf^^i^(kinr.^,.)} 



(9) 



as the normalisation constant. A peculiarity arises with this 
choice of boundary condition; it can be shown (Fisher et al. 
1994c) that the mean value of a field with the expansion (4) 
must be zero. This is dealt with by adding the mean value 
to the sum A (^A) -\- X^;^^. Since we shall be concerned 
with fluctuations around the mean, the extra term cancels. 



3.4 Compact Notation 

In later equations, the explicit reference to Imn when refer- 
ring to a particular mode is somewhat cumbersome. Hence, 
we shall on occasion use Greek indices as a shorthand for 
modes, keeping in mind the implicit relation /t {I,m,n). 
This is related to other representations of scalar quantities 
by 

Aimn = {A\lmn) = A^. (10) 



4 REDSHIFT SPACE DISTORTIONS 

4.1 Transformations to redshift space for 
inhomogeneous samples 

Having laid out the basic transformations, we shall now turn 
our attention to the relationship between the transform of 
the observed density field, and the true density field. The 
former is modified in two ways: firstly, because of the selec- 
tion function, the mean density of galaxies drops off with 
distance; secondly, the observed distances are redshift dis- 
tances, not real distances. A further complication arises if we 
wish to weight the galaxies depending on their redshift (in 
some applications, inverse selection function weighting may 
be appropriate, e.g. in reconstructions of the real-space den- 
sity field. Fisher et al. 1994c). We will identify the observed 
coefficients by pimn, and the coefficients of the underlying 
real density field by pimn- 

In general the observed redshift density field can be 
defined as a sum of Dirac delta functions. 



P(s) = ^fc(s(r, v) 



)fc(v 



(11) 



projected from the phase space distribution of galaxies. This 
can also be expressed in terms of fluctuations about a back- 
ground mean density field, which will in general be spatially 
varying due to selection criteria such as flux limits and un- 
observed regions of sky. We shall denote this generalised 
background by /9o(r). Hence the observable density field is 
now expressible as 



p(s) = po(r)(l + «(s)). 



(12) 



The background density, po , may in fact be a function of the 
real and redshift-space coordinate systems, as the condition 
for inclusion into a catalogue may be functions of either sys- 
tem. Here we shall only consider real-space selection criteria. 
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Applying the harmonic transform to the redshift density 
field, we construct 



Cln / w{s)p{s,)jl{kl„s)Ylm{d,4>)d S 



(13) 



where we have included an optional weighting function w(s). 
Using the continuity condition (f" s p(s) = (f'r p(v), we can 
transform the integral to real space; 



Plmn = Cln / w{s)p{r)ji{kins)Yim{9,(l>)d^r. 



(14) 



where the distortion is now restricted to the argument of 
ji{kins) and the weighting function w{s). We see also that 
our choice of boundary condition has avoided the introduc- 
tion of surface terms into this expression. In general, we 
should expect that these would be negligible in any case, 
being dominated by the effects of galaxy selection near the 
survey boundary. 

While expression (14) is exact it does not represent a 
complete solution to our problem, as the relationship be- 
tween s and r has not been supplied. To complete the anal- 
ysis we shall proceed by perturbation expansion oi ji{ki„s) 
to first order. From Section 2 we see that the velocity field 
can be related to the gravitational potential by a gradient 
operator, and that the redshift displacement is the radial 
projection of this. Using our results (Appendix A) for the 
vector spherical harmonics we see that 



-r • V 



-«p(r: 



^o'^ ^ Clnki^ 8lm„j'l(klnr)Ylm(d , (f)) , (15) 



is the radial potential field term of the velocity field and 
where ji(z) = dji/dz. Simn is the transform of the fractional 
galaxy number overdensity field 8{r), and we shall assume 
that the galaxies are biased tracers of the density field, in the 
sense that 8imn = b{8imn)p- This condition is only assumed 
to apply for the modes used in the analysis. In particular, no 
assumptions need to be made about the short-wavelength, 
highly nonlinear modes, as these will not significantly af- 
fect the long-wavelength modes analysed here. With this 
assumption, Qq^ ^ /? = Qg^/fc, if we assume that b is inde- 
pendent of scale k. For I > 0, 8imn is simply related to the 
transform of the density field pimn by pimn = p8imn, where 
p is the mean number density''' . 

Note that the modes are independent in real space: 



{8lmn8*i m' n') — P (kln)8fii 8^^i 8^„i 



(16) 



where P{k) is the power spectrum of mass fluctuations, and 
8^ is the Kronecker delta function. 

Returning to the expansion of the observed redshift 
space density field, we can expand ji{ki„s) to first order 



by 

]l(klriS) 



Ji(ki„r) + u(v)-^ji(kinT) 



+ Interestingly, we do not need to assume gravitational instabil- 
ity; we simply need to determine the velocity field for which the 
real-space density power spectrum is isotropic. This is discussed 
in Section 6. 



— Jl{kinr) + 

c;/„/ ^'j"^'"' j'v(ki,„,r)j'i(kinr)Yii,n' ■ 

I' m' n' 

(17) 

Consider now the effects of a inhomogeneous background 
density. Expanding out the scalar perturbations, 8{v), and 
using the orthonormality oi ji{ki„r) and Yim{d,(t') (see Ap- 
pendix C), we find that 

Plmn = (po)lmn + ^ i^Wnn' + I^Vi'^nn' )^l'm'n' , (18) 
I' m' n' 

where 

Imn — Cln j pO {r}w{r}ji{kinr}Yim{e,<t>)dr (19) 

is the transform of the observational, weighted mean den- 
sity field. In the general expression for pimn the transition 
matrices are defined as 

^wZn' = {lmn\pow\l'm'n'), 

= cinCVn' j po(r)w)(r) X 

Jl(klnT)Ylm(eA)jv(h'n'r)Yirn,'(e,^)d\, (20) 

and 



Trmm 
^IV nn' 



y Po(r)^[«;(r)j,(A;,„r)] 

y.Yim{e, <t>)j'i,{ki,n'r}Yifn,'{0, 4>)d\. (21) 

In a more compact notation, we arrive at a final expres- 
sion relating the undisturbed density field with the observ- 
able redshift distorted density perturbations, D = {pimn — 

{po)lmn)/p; 

Di. = Y,K^-' (22) 
where 

a;: = $;: + /?t/;. (23) 

Equation (22), in this its most general form, is one of the 
main results of this paper. In comparison with the Kaiser 
equation, we note that redshift distortions causes the mix- 
ing of modes, as we expected to happen from the general 
grounds discussed in the introduction. This expression is ex- 
act to first order, and correctly deals with the finite extent 
of the redshift survey. 

4.2 Division into Radial and Angular Selection 

One of the main advantages of using spherical coordinates is 
the natural division into radially dependent and angularly 
dependent functions. In application to large angle redshift 
surveys this is particularly useful. Thus we can decompose 
/9o(r) into 



po{r) = W{0,<f,)po{r) 



(24) 



where W{9,(l>) is an angular selection function, and /9o(j') 
now embodies the radial selection function. 

In the case of all-sky redshift surveys such as the IRAS 
QDOT and Berkeley 2-Jy and 1.2-Jy surveys, there is a re- 
gion around the galactic plane, the zone of avoidance, that 



6 A.F. Heavens and A.N. Taylor 



is excluded from the survey due to source confusion. This 
masking introduces, for example, a large quadrupole term in 
to the distribution of sources, which is then mixed with other 
modes by the redshift-space distortion. With this separable 
form for the selection function, the transition matrices be- 
come (cf Peebles 1973, Scharf et al. 1992) 



wr, 



T T-nn' 



(25) 



(the summation convention is suspended) and 

cinCVn' I Po{r)w{r)jii{kii„ir)ji{ki„r)r'^dr 



V"', 

' nn' 



cincvni I pQ(r)-^ [w)(r)ji/(A;i/„/r)] j,'(A;i„r)r^(ir 
ar 

Wl^r = I Yvm' M{U)Yil,dU, ^ = {e,4>). (26) 



For a complete all-sky survey, the expressions simplify, since 
W reduces to delta-functions: PF,'^™ = 8^,8^^,. So, in the 
absence of a mask, the mixing is between modes with the 
same I and m, with considerable simplification. 



4.3 The Covariance Matrix of the Dimn''s 

In the general case, we can compare theoretical models with 
the observations by constructing ensemble averages; 



5D 



(D,.D.) JK + HKVK + ?V:tP{i.t + A"„ (27) 



where, in practice, we split the coefficients into real and 
imaginary parts, which are independent provided that the 
mask mixing matrices W are real (this is true for a zone of 
avoidance mask, and very nearly true for the IRAS mask). 
This explains the factor |- in the power. The last term is the 
shot noise, which can be obtained by the methods of Peebles 
(1973); 



A, 



cinCii„' / po{r)w^ {r)ji{kinr)n' {ki'„' r)r^ dr x 



T^(Y,min))M(n)T^(Y,f^,(n))dn. (28) 



Here 'Pjn.i/ represent real or imaginary parts, depending 
on whether -Dju.i/ are real or imaginary. In order to cal- 
culate the shot noise contribution, we use the fact that 
Tesseral harmonics (spherical harmonics of the type TZe or 
Im[i^'"(Q)]) may be written, for example, 7^e[i^'"(Q)] = 
|'[Yim(f^) + Yim(^)]) which allows us to calculate the angu- 
lar integral in terms of four W^;™™ matrices. 

Equation (27) will form the basis of estimation of pa- 
rameters via likelihood methods, but it is instructive to 
pause and consider it for a moment. Essentially we make the 
comparison in the observational domain, assuming a power 
spectrum P{k) and a value of /?. Notice that the observed 
quantities depend on both of these, so there is hope that 
a maximum likelihood method would yield information on 
the power spectrum and on /?. Indeed, for certain weight- 
ing functions, the mixing matrices are not very broad, so 
estimation of the power spectrum may be possible without 
great uncertainty. This will be the subject of a subsequent 



paper; here we assume a parameterised form for the power 
spectrum. 



4.4 Small— Scale Velocity Structure 

The small-scale velocity field of galaxies introduces an ex- 
tra fine structure in redshift space. Once overdensities col- 
lapse and subsequently virialize, their velocity distribution 
relaxes towards a Maxwellian distribution. Galaxies in these 
structures are scattered radially about their true positions in 
the cluster producing extended structures in redshift-space. 
These have been called "Fingers of God" due to their radial 
alignment. This random behaviour is strongly nonlinear in 
origin, but works in the opposite sense to the linear com- 
pression of large-scale features. Small-scale structures are 
smeared out, resulting in the strong aliasing of power from 
the higher to smaller, linear wavenumbers. 

We can mimic the effect of the small-scale velocity struc- 
ture in redshift-space by randomly assigning each galaxy an 
additional peculiar velocity drawn from a Maxwellian distri- 
bution (cf Peacock & Dodds 1994), and projecting radially. 
This is approximate, since it assumes that the velocity dis- 
persion is independent of location. Formally the effect can be 
calculated by first introducing a small, random perturbation 
into the position of each galaxy, and averaging: 

i)'(r) = (po(r)«[s'(r)])„ (29) 

where the redshift-space coordinates are displaced by 

s(r) — r s'(r) = s(r) -|- e(r). (30) 

Here e(r) is a Maxwellian distributed random variable, with 
(e(r)) = and (e,(ri)ej(r2)) = (T^5,j5_D(ri - r2), where aH 
is the assumed uncorrelated 3-D velocity dispersion. 

It is straightforward to show that the modes of the per- 
turbed observable galaxy density are related to the unper- 
turbed field in redshift-space by a convolution. We shall 
leave the outline of this to Appendix B, but quote the re- 
sult, which is 

Dlmn = nn'Dlim'n' , (31) 

where the cluster scattering matrix is given by 



mm 



T/2 



P(l" - y)jl{klnr)Ylm{Or, 4>r 



'XJv{h'n'y)Yvm'{Sy, 4>y)r^dr dUr y^dydUy , 



(32) 

To find the scattering matrix suitable for projections along 
the radial direction, we decompose the Maxwellian velocity 
distribution, p(r — y), into a radial and angular dependent 
part; 



p(r ■ 



-(r-y)2/(2<T2) 

(277)3/2^3 
„-(r-j,)=/(2<T=) 



27r(7 



-ry(l-l-i')/a 



27r(72 fry 





' 1 " 




_ry_ 



(33) 



where we identify the angular dependent term in the ex- 
pansion as a negative exponential distribution. Allowing the 
width of this term to vary as cr ^ cr', we can collapse this 
distribution into a pencil beam. As a' 0, we find 
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-ry(l-ii)l<y' 

l^Ku'^ Iry 



fe(l-/*)/(27r). 



(34) 



Applying this constraint to the scattering matrix equation 
(33), we arrive at the symmetric matrix for a radially pro- 
jected random perturbation in galaxy positions; 



cK cK 



27r(7 



-Ji(hnT)]v(kvn'y)rdT ydy. (35) 



Note that we absorb the mask into the translation from 
real to (unsmoothed) redshift space, so the smoothing is 
a delta function in m and m'. This random scattering of the 
observed density field can be folded into our calculations for 
the real underlying modes of the density field by the final 
expression 



(36) 



where summation is implied. 

Having specialised to the case of radial random fluc- 
tuations in the positions of galaxies, it apparent that the 
effects of measurement errors in galaxy redshifts can also 
be accounted for by adding this to the velocity variance in 
quadrature. 

4.5 Low-order modes and uncertainties in the 
Local Group Velocity 

The relationship between the overdensity coefficients 8imn 
and the density coefficients is simply pimn = pSimn, except 
for the I = modes, for which the radial selection function is 
involved (see Fisher et al. f 994c for more discussion of this 
point). In addition, if we change the frame in which we mea- 
sure the peculiar velocities (perhaps by choosing the Local 
Group frame, and revising the Local Group velocity with re- 
spect to the microwave background), the coefficients remain 
unchanged, except for the dipole term. For these reasons, 
it is sensible to restrict one's attention to the quadrupole 
term and higher. If the survey is not all-sky, changing frames 
does affect I ^ 1 modes, since there is some mixing between 
modes, but these are small if the masked regions of the sur- 
vey are small. 



5 STOCHASTIC PROPERTIES 

Having set out the formal relationship between individual 
harmonic modes of the underlying density field in real-space, 
and that of a set of modes measured in an observationally- 
limited redshift survey, we now move on to discuss the statis- 
tical properties of the observable density modes. Our main 
result here is that the observable density field no longer re- 
tains all the statistical properties of the true density field. In 
particular the usual assumption of statistical homogeneity 
and isotropy of the density field break down for a finite red- 
shift survey. It is this fact that allows us to use a statistical 
approach to determine the amplitude of the distortion, and 
so measure /?. 

In Section 4.3, we calculated the covariance matrix of 
perturbations in redshift space, and related it to the power 
spectrum of the true density field. Using this relation, we 



may make a random phase hypothesis regarding the phases 
of the true density perturbations. This leads us to regard 
the harmonic modes of the true density field as Gaussian 
distributed. A fundamental property of Gaussian random 
fields is that the linear transform of the density field still 
obeys the Gaussian form. In Section 6, we shall use this 
property to measure the degree of distortion from statisti- 
cal isotropy and homogeneity present in each mode of the 
observable density field. 



5.1 Statistics of Observable density fields 

The lowest order non-trivial statistic of the perturbations 
is the covariance matrix of each mode of the density field, 
{DimnDii^i„i), where angled brackets, (), denote ensemble 
averaging, and is given by equation (27). The last term in 
that equation is due to discrete shot noise. Clearly, as we 
shall assume that this is a random sampling of the density 
field by galaxies, this term is independent of redshift space 
distortion. However it is still affected by the finite volume 
and both angular and radial selection functions. 



5.1.1 Recovery of Kaiser's result for k ^ 

It is interesting to consider the relationship between the av- 
erage redshift space power P''{k) and the real-space power 
P{k), and to compare it with the analogous result of Kaiser 
(1987) for the distant-observer approximation {P''{k) = 
[1 + 2/?/3 -I- /]^/5]P{kj). Since the redshift distortion in- 
troduces mixing between modes, the ratio of P''{k) to P{k) 
depends on the power spectrum, in contrast to Kaiser's anal- 
ysis: 



P'ik) = {\D,f) 



m + i^v,yp{k^) 



(37) 



where the averaging is done over bins in wavenumber, and 
we ignore shot noise, incomplete sky coverage and incoher- 
ent velocity dispersion. Multiplying out the square, and per- 
forming the summation, for a flat power spectrum, {P{k) = 
constant), we calculate the coefficients of /? and /9^. These 
are plotted in Figure 1, along with the values 2/3 and 1/5 
from Kaiser's analysis. 

Two effects on the values of the coefficients that are 
not seen in the Kaiser terms become apparent. Firstly, it 
is apparent that edge effects at low ki„ attenuate the am- 
plitude of the modes. As we only have a finite sample, this 
expresses the fact that information about modes larger than 
the survey are not directly measurable. 

The second effect is due to the small-angle approxima- 
tion used to derive the Kaiser equation. This has the effect of 
pushing the range of validity of the Kaiser equation to higher 
kin and I modes. This is the same as saying that a mode in 
the linear regime subtending a small angle will have to be at 
a correspondingly larger distance. Both of these effect con- 
spire to push the range of validity of the Kaiser equation 
to smaller scales, where we expect nonlinear coupling in the 
growth of structure to become dominant. If this is the case, 
then the range of applicability of the Kaiser equation may 
be somewhat limited. 
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5.2 Likelihood Analysis of the Observable density 
field 

A natural condition for density perturbations formed in tlie 
early universe to satisfy is that of random phases. The Cen- 
tral Limit Theorem implies that these density fluctuations 
will be randomly Gaussian distributed, with all the neces- 
sary statistical information contained in the covariance, or 
correlation, function. In the case of a continuum density 
field, which has an infinite number of degrees of freedom, 
the distribution function is dependent on the values of the 
density field at all points, and is called a distribution func- 
tional. 

The distribution functional for density perturbations 
expresses the probability for a given configuration of fluctu- 
ations to occur. This configuration can be quantified either 
in terms of the field amplitude at each point, or more com- 
pactly, in terms of the amplitude of each harmonic mode. 
Generalizing the case of a multivariate Gaussian, the prob- 
ability distribution functional of the density field can be ex- 
pressed as 



P[8] = {detA)~^'^ exp ^-i J 8{r)A~^ (r, y)8{y) d\d' 



(38) 



phere 



(detAY'^ = J D[S]exp ^-^ J S(r)A~^(r, y)S(y) d^d^ 

(39) 

is the natural extension of the determinant to an infinite- 
dimensional function space. The Wiener measure, i)[5(x)], 
is defined in the limit as limjv^oo fl^-i d8{yin) I {^-k)""^^ . The 
function A(r,y) = {(|r — y\) is the correlation function of 
the density perturbations. An equivalent expression can be 
found by expanding the density field in terms of its harmon- 
ics, leading to the expression 



P[8] = (detyl)"^/^exp 



/ _ 1 8lmn8*im 

\ 2^ P(k,n) 

\ Imn 



(40) 



where we have used the statistical independence of the har- 
monic modes to simplify the summation. 

The amplitudes |5imn| follow a Rayleigh distribution, 
and the phases are random, or equivalently, the real and 
imaginary parts of Simn are independently Gaussian dis- 
tributed. We can use this to property to split the probability 
distribution functional of the density field into two indepen- 
dent parts, one for real terms and one for the imaginary 
terms; 



P[8] = P[7^e(^^„„)] P[lm{8im„)l 



P[7^e(^^„„)] = (detyl,je)"'^'exp 



P[Im{8imn}] = {detAxm} exp 




E Im{8lmn) 
2P(k,n) 

Imn 



A general property of Gaussian random fields is that 
linear operations on the random variable are also Gaussian 
distributed. Our transformation equation (22) from the real 
density to the observed density field is of a linear form, and 
so the real and imaginary parts of Dimn are also Gaussian 
distributed. The linearity of equation (22) is essential here, 
and arises as a results of choosing spherical harmonics as 
our basis for expanding perturbations. In other bases, both 
the redshift projection and window function may be nonlin- 
ear operations, and our assumption of Gaussianity would be 
doubtful. 

With this in mind, we are now in a position to construct 
a likelihood functional for the observed density field. Apply- 
ing the transform (22) to the underlying density field, and 
using equation (27) for the correlation matrix of observable 
density modes, we find 

C[D\l],P{kin}] = (detC,je)"'^'(detCi„)-'/' 

X exp I -^'^Tle{Dimn){Cl^''^„,j^^)~^Tle{Di,m'n') 1 

\ Imn / 




mm' —1 
Tfl{ D lmn){C III -j^^ 



where Atzc and Azm are the correlation matrices of the real 
and imaginary coefficients respectively. 



X exp 



(42) 

as obtained by Scharf & Lahav f 993 and Fisher et al. f 994a. 
The inhomogeneity of the observed density field is manifest 
in the likelihood by the correlations between modes, which 
also contain valuable phase information. 

In principle we can proceed in a number of ways; as- 
sume some functional form for the true power spectra, and 
maximize with respect to /?, or allow some freedom in the 
shape of the spectrum and maximize C for a fixed /?, or al- 
low both the shape and /? to be determined simultaneously. 
In the next section we discuss and apply some of these pos- 
sibilities. 

It is worth noting at this point that the method we 
have outlined is in essence a power spectrum analysis. Aside 
from allowing us to represent compactly the amplitude and 
phase information of the density field in terms of fundamen- 
tal modes, it allows us to separate the density perturbations 
into linear and nonlinear regimes. This is an essential fea- 
ture, as the method outlined is exact only to first order. Even 
the inclusion of a correction term from the random velocities 
found in clusters only serves to suppress the contamination 
of small scale modes aliasing as linear modes. This high- 
lights the general problem of nonlinear modes coupling with 
linear modes, due to selection functions mixed with redshift 
distortions. With careful choice of weighting functions for 
the data, the mixing can be made minimal. This contrasts 
strongly with methods based on the extraction of redshift 
space distortion from galaxy correlations, where the distinc- 
tion between linear and nonlinear regimes is ambiguous. 

However, a general feature of the statistical approach 
is the degeneracy that will arise in the parameter space of 
redshift distortion, /?, and the intrinsic amplitude of density 
distortions. An increase in the distortion has the same effect 
on radial modes as that of increasing the amplitude of ra- 
dial perturbations. However, homogeneity demands that an 
increase in power is distributed equally between all modes, 
regardless of phase. Hence the likelihood functional, (42), 
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can, to some extent distinguish between the two cases, and 
break this degeneracy. 



5.3 Optimal Weighting of the data 

In Section 4.2 we allowed ourselves the option of an ar- 
bitrary weighting scheme for each galaxy in the sample. 
This is a desirable feature, as is it not apparent that each 
galaxy in a sample subject to selection and bias contains an 
equal amount of information about the underlying density 
field. Indeed, it turns out to be essential to adopt a suitable 
weighting scheme, to avoid the convolutions arising from the 
selection function introducing coupling between too many 
different modes, which becomes computationally difficult to 
handle. 

Given that we do have prior information about how 
the sample was selected, from the magnitude limit, window 
function, etc, it is reasonable that we use this information 
to optimally weight each galaxy when measuring a specific 
statistic. 

A weighting scheme can be called optimal if it min- 
imises the variance in the measured statistic, in which case 
it is called a minimum variance weighting scheme. Such a 
minimum in the variance of a statistic can be found by solv- 
ing the functional equation 



8w{y 



■Var[/(r)] = 0, 



(43) 



where / is the statistic being measured. The variance of f{r) 
can be calculated directly from the likelihood function using 
the approximate relation 



Var[/(r)] 



6' In C[f{y)] 
8f{rf 



(44) 



where — In £ is the error functional and — f/Var[/] is called 
the amount of information. Applying this expression to the 
Gaussian random field model, we find 



Var[/(r)] ~ 2 ( 



(45) 



where C is the covariance matrix and C" = hCjhf is the 
functional derivative of C with respect to /(r). 

For the moment we shall ignore the effects of redshift- 
space distortion, and concentrate on finding the minimum 
variance weighting function for a function that is distortion 
independent. After this we shall be in a better position to 
generalize the problem to include distortions. 

Taking the diagonal elements of the covariance matrix 
only, Cii = diagCjj, and assuming that the power spectrum 
is approximately constant over the range of interest, we find 
that 



w\r)^\r) [ — + P{ki„) ] jfiki„r) d\, (46) 



where we have used the closure relation for spherical Bessel 
functions (Appendix C) to reduce the summation over 
modes. 

The two statistics we are most interested in obtaining 
from redshift space are the power spectrum of density fluctu- 
ations, P{ki„), and the distortion parameter, /?. In the limit 
of no distortions, only the power spectrum is relevant. The 



variance of the power spectrum, measured over some small 
range in harmonics, is 

VarP(A;,„) 



2 / w^{r}ct,^{r} 



1 + 



jt(kinr) d'-'r 



I J w^{r)4>^{r)j^{ki„r) d^r\ 



(47) 



Taking the functional derivative with respect to the weight- 
ing function, w{y), and remembering the functional relation 
8w{x) / 8w{y) = 5_d(x — y), the minimum variance condition 
is 



■In 



VarP 



8w P2 

4 / w^r)<f>^{r) 



1 + 



jt(kinr) d'-'r 



^ 4 / w(r)<f(r)Pi(ki„r) d\ 
J w^(r)cl>^(r)jf(kinr}d-^r 

with solution 

f 



3f{ki„r) d^r 



w[r) 



1 + P{kin)<f>{r) 



(48) 



(49) 



where we have ignored an unimportant amplitude. This co- 
incides with the optimal weighting function derived some- 
what differently by Feldman, Kaiser & Peacock (1993), al- 
though both derivations make use of the properties of the 
Gaussianity of the density field. This weighting scheme has 
the expected form of signal/ (signal -|- noise) that one would 
expect from discrete data that can be decomposed into a 
Gaussian distributed signal and Poissonian shot noise. A 
big advantage of this weighting function is that the mixing 
matrices ^^'^n' much narrower in k-space than those 
for equal weighting, and fewer modes need to be considered 
(see figure 2). This form leads us to expect that a mini- 
mum variance weighting scheme can be found by applying 
the approximation 



Var[/(r)] 
Pir) 



dlnC 



dlnf 



(50) 



where C is approximated by equation (46). This allows us 
to extend the analysis to the case where P ^ 0. The prob- 
lem is now much more difficult, since the /? term involves 
derivatives of the Bessel functions. In order to proceed, we 
assume that these derivatives obey an approximate orthog- 
onality relation; since this is not quite true, the weighting 
functions are not quite optimal, but may still be useful. We 
assume further that w{r) is locally a smooth power law with 
slope 7: w{r) ~ (r/ro)"'. We can then repeat the derivation 
of of the weighting function for the power spectrum, to find 



wp{r) 



{1 + 1^ + yl^ro/r) 



(51) 



1 + (1 + P + iPro/ry P{kin)(t>{r) 

which in the special case of /? = reduces to equation (49). 
Applying the same formalism to the variance on /? yields 



Wi3{r 



+ P + jPro/r){l + jro/r) 



(52) 



1 + (1 + 13 + -/iSro/ry P(kin)(t>(r)' 
In Figure 3(a) we plot the selection functions for QDOT, 
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the 1.2-Jy and 2-Jy IRAS surveys. Assuming that the power 
spectrum of perturbations can be approximated by a mean 
value, P{k), we show in Figure 3(b) the weighting functions 
corresponding to the 1.2-Jy catalogue. The solid curve shows 
uniform weighting, the Feldman, Kaiser, Peacock weighting 
scheme is shown by a dashed curve, while dot-dashed and 
dotted curves show the redshift-space corrected weighting 
schemes. All curves are normalized to unity at r = 3 x 10* 
km~^s. 

Interestingly, the optimal weighting schemes for P{k) 
and P are very similar, and significantly different from the 
uncorrected weights. Two effects are can be seem. At low 
redshift, the jflro/r terms are dominant, and correct for 
the misplacing of the weights in redshift space. Beyond To, 
the dominant effect is the amplification of the power by a 
factor 1 + /?. 

It is worthwhile clarifying some of the approximations 
used in deriving the two weighting functions, wp{r), and 
Wf:i{r), and what they actually correct for. Our main assump- 
tion is that the minimum variance weights will be approx- 
imately of the form ~ (C'/C)^, which is true for Gaussian 
distributed variables. Also we have assumed orthogonality 
between the derivatives of the spherical Bessel functions in 
order to derive equations (51) and (52). While this allows 
us to introduce factors oi 1 + P + jflro/r, we have ignored 
the angular dependence of the distortion. However, as our 
weighting scheme is arbitrary, we are free to weight galax- 
ies as we please, and so are only concerned with finding an 
approximately optimal weighting scheme. Its value can then 
be assessed by its influence on the measured statistics. 

Given these caveats, the weighting schemes, equations 
(51) and (52), optimally weight galaxies that randomly sam- 
ple an underlying density distribution, subject to both a se- 
lection criteria in real space and a distortion that shifts the 
radial coordinate of each galaxy. The weights correct statis- 
tically for the distortion, selection and shift when measuring 
a specific statistic of the density field. 

We shall return to the effects of these these weighting 
schemes on the transition matrices later. 



6 APPLICATIONS TO SIMULATED AND 
REAL DATA 

6.1 N-body simulations 

We have tested the method using a random subsample of an 
N-body simulation provided by Bertschinger & Gelb (1991 
for details). We choose an observer at random, and form 
the observer's redshift space map within a large sphere. 
In order to make the sphere radius large compared to the 
scale at which the fluctuations become nonlinear, we repeat 
the box before constructing the redshift space map. In the 
simulation we used, the power spectrum in real space is 
rather accurately given by P(k) = [k / koY '^'^ / k^) , with 
ko = 0.45/t Mpc~ . In the first test, we analyse the real-space 
map, constructing likelihood contours for the quantities /? 
and ko. In essence the real-space map has /? = 0. The re- 
sults are shown in Figure 4, which shows the maximum like- 
lihood solution to be an unbiased measure of both ko and /?. 
In Figure 5 we analyse the redshift space map, constraining 
the power spectrum to have the right form. The maximum 
likelihood solution lies close to the expected /? = 1. 



6.2 IRAS 1.2 Jy sample 

We have applies the method to the IRAS 1.2 Jy sam- 
ple, containing ~ 5000 galaxies. We form a catalogue with 

' max 

= 200/t"^ Mpc and a conservative galactic plane cut 
of \b\ > 10°, to minimise any residual confusion with galac- 
tic sources. This leaves 4511 galaxies, with a mean observed 
density (Fisher, private communication) of 



Po 



ar-^'^{l + {r/r,f)-^' 



(53) 



with a = 0.189, a = 0.421, /?, = 1.913 and r, = 50.1. Units 
are Mpc. We assume a parameterised form for the power 
spectrum given by Peacock & Dodds (1994): 



m 



27rM(fc/fco) 

P [1 + (k/kc) 



(54) 



where we take the amplitude yl to be a free parameter, ko = 
0.29, kc = 0.039/t Mpc~^. A is related to the variance in 
the real-space number overdensity in spheres of radius 8h~^ 
Mpc by as, IRAS = 0.69 \/A. (Note that as is often a linear 
theory extrapolation; the value here is calculated from the 
assumed present-day power spectrum). 

We treat /? and asjRAS as free parameters, and present 
results including modes up to A; = 0.1/tMpc"\ This gives 
a total of 604 modes, with maximum values I = 17 and 
ra = 6. The convolutions contain modes up to I = 30 and 
n = 20. Shot noise is not very significant for these modes, 
with signal-to-noise rarely less than, and usually far in excess 
of, ten, even for the highest I modes considered. 

Likelihood contours are shown in the asjRAS — P plane, 
in Figure 6 for the weighting scheme designed to minimise 
variance in the measured power. It is impractical to choose a 
different weighting scheme for each wavenumber, so we take 
a single power equal to 6000 h~'^ Mpc'^, roughly what we 
expect at a wavenumber of 0.15 h Mpc~^. The maximum 
likelihood values are /? = 1.1 and asjRAS = 0.68, in agree- 
ment with the results of Fisher et al. 1994a. Note that this 
graph also contains the redshift space r.m.s., from the limit 
/? ^ 0, as, IRAS, z = 0.87 ± 0.05, as found by Feldman et al. 
1994 for the QDOT sample, and marginally higher than the 
value of 0.80 found by Fisher et al. 1993 for the 1.2-Jy sam- 
ple. The error analysis is complicated because the likelihood 
contours are very elongated. The ridge of likelihood extend- 
ing towards high /? arises from degeneracy between /? and 
the amplitude of the power spectrum, which arises if the red- 
shift space power is much larger than the real-space power 
because of large distortion, as is clear from examination of 
(27) in the high-/? limit. 

The formal error from 8ln C = —0.5 projected on to the 
parameter axes gives the formal errors: 



13 = 1.1 ±0.3 
as,iRAS = 0.68 ±0.05. 



(55) 
(56) 



The errors are perhaps unduly optimistic. The ridge of like- 
lihood extending towards as,iRAS = 0.87 does not exclude 
/? = with confidence, the 95% confidence level being at 
8lnC = —3 for the Gaussian approximation (Press et al. 
1992). If we were to increase the wavenumber range so that 
nonlinear effects became important, we would expect the 
maximum likelihood value of /? to decrease. We begin to see 
this when we increase kmax to 0.12; the maximum likelihood 
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value drops to 0.5, and the results become considerably nois- 
ier than those in Fig. 6. 

One of the powerful aspects of this technique is that 
we have control over the wavenumbers involved, and this 
wavenumber limit allows us to avoid nonlinearities. Peculiar 
velocities of the order of 300 km s~^ will affect wavenumbers 
down to ~ 0.3, which is also roughly where the density field 
becomes nonlinear. For this reason we expect the correction 
for wavenumbers less than O.f to be small, and this is indeed 
the case. Figure 6(b) includes the correction, 6(a) does not. 



6.3 The value of Qo, with COBE, inflation and the 
abundance of rich clusters 

An interesting results is obtained if we use measurements 
of the microwave background made with the Cosmic Back- 
ground Explorer (COBE) satellite (Smoot et al. f992), and 
make the further assumptions, motivated by inflation, that 
the power spectrum P{k) oc k sls k ^ 0, and that the Uni- 
verse is spatially flat, which requires a cosmological constant 
if Qo < 1. The QjJ depend ence of the velocities is virtually 
unchanged by the presence of a cosmological constant (La- 
hav et al. f 99f ) and gravitational instability theory allows us 
to relate the amplitude of the power spectrum to the COBE 
temperature fluctuations (Efstathiou et al. f992): 



Pik)^Bk;B = ^^^' 



5 V T 



2C_Y „-1^54 f h 

Ho) ° Impc 



(57) 



With a quadrupole of f7.0 x fO"'' K (Banday et al. f994) 
and T = 2.735K, we get 



P{k) = 5.93 X fO^fcQ^^^^" 



f.7 X fO-5, 

with k measured in h Mpc~^ . Relating this to our parame- 
terised power spectrum, we find the bias parameter, assumed 
independent of scale, at least on these large scales: 



b = 0.842Vl'no " 



f. 22(78, lilAsf^o'''' 



(59) 



We can also incorporate the constraint found by White et 
al. (f 993) from the abundance of rich clusters: 



(78, p ~ 0.57 Qq 



(60) 



with a quoted error of about fO%. We can combine these 
data to provide weak constraints on the parameters Qo and 
b. This will be the subject of a forthcoming paper (Taylor 
& Heavens, in preparation), but preliminary results are that 
Qo >0.2 and b > 0.35. 



6.4 Peculiar velocity power spectrum 

Although we have assumed gravitational instability for most 
of this analysis, it is not strictly necessary. Under less re- 
strictive assumptions, we can calculate the peculiar velocity 
power spectrum required to make the observed density data 
consistent with an underlying isotropic field. We need to 
make the following assumptions explicitly: the velocity field 
is curl-free; galaxies obey the continuity equation; modes 
grow independently of k. 

In this case, we may write the peculiar velocity in terms 
of a velocity potential, v(r) = — V^*, and the linear con- 
tinuity equation in expanding coordinates gives, V^^* = 



a{t)d8 1 dt. Assuming modes grow independently of k gives 
a{t)d8/dt = F{t)8. The power spectrum for the velocity 
is then (|vj,p) = F (ty (\Sy^\^ \) / k^ . For gravitational insta- 
bility, F(t) is just 13/ Ho, so by determining jS and ) 
by likelihood techniques we determine the velocity power 
spectrum, independently of the assumption of gravitational 
instability: 



too- 



i{k/ko) 



k^ [1 + (k/kc)-^ ''] 



(kms 



(61) 



defined such that the variance in the three-dimensional pe- 
culiar velocity field is (7^ = J A^dlnk. ko = 0.29 and 
kc = 0.039/t Mpc~^ are obtained from the assumed power 
spectrum shape. Since (yl/?^)max — 1, we obtain the required 
velocity power spectrum, shown in Figure 7. Note, that, al- 
though the expression (6f) involves /?, the calculated is 
independent of the assumption of gravitational instability. 



7 DISCUSSION AND FUTURE PROSPECTS 

We have presented here a very powerful method for extract- 
ing the quantity /? = Qg^/fc and the power spectrum from 
near all-sky galaxy catalogues with the minimum of assump- 
tions. In essence, the degree of anisotropy in the redshift 
space catalogue is exploited to measure the distortion from 
the expected real-space isotropy. The method is exact to 
first-order, so represents a significant advance over Fourier 
methods. In addition, as a power spectrum approach, we can 
limit our application easily to scales in which the linear ap- 
proximation holds, one of the reasons why it is superior to 
correlation function studies. The method exploits the pres- 
ence of peculiar velocities, but does not rely on actual mea- 
surements of peculiar velocities, and thus is not sensitive to 
any uncertainties over how well they can be measured. We 
believe that this formalism represents the best method for 
analysing near all-sky galaxy redshift catalogues. 

It is worth comparing the methods presented here and 
in Fisher et al. f994a. In the latter case, the galaxies are 
weighted with arbitrary functions of redshift, and projected 
on to the sky. The resulting sky distributions are then anal- 
ysed with spherical harmonics. This has an advantage over 
the method of this paper in that one can try to tailor the 
weighting functions to maximise the effect of redshift distor- 
tions. In one sense, the work in this paper can be regarded as 
a particular choice for the set of weighting functions, but this 
misses the power of the method. The expansion in spherical 
harmonics is analogous to a full Fourier expansion in carte- 
sian coordinates, offering similar advantages such as direct 
measurement of the power spectrum (albeit convolved with 
some window functions) and direct analysis of anisotropy 
of the coefficients. An obvious extension to the application 
discussed here is to allow the form of the power spectrum 
to be free. This work is in progress, and should allow a si- 
multaneous determination of P{k) and /?. A further appli- 
cation of this technique is to surveys which are not all-sky. 
Without all-sky coverage, there is mixing of different modes 
with different I and m. This makes the computation of like- 
lihoods computationally expensive, but the complications 
can be minimised if the survey has azimuthal symmetry. In 
this case mixing only occurs between modes of the same m 
value, and the resulting likelihood matrices become block- 
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diagonal. This is an important simplification which makes 
the computations tractable. 

Our analysis of the IRAS 1.2 Jy galaxy catalogue yields 
the following parameters: /? ~ 1.1±0.3, and a real-space am- 
plitude CTgjRAS — 0.68 ± 0.05, confirming values found by 
other methods (e.g. Fisher et al. 1994a, Peacock & Dodds 
1994, Dekel et al. 1993). However, with the exception of 
Fisher et al. 1994a, the other methods which do not use the 
peculiar velocity field directly have some further approxi- 
mations which are hard to justify, or use correlation func- 
tion distortions, which are hard to separate into linear and 
nonlinear regimes, and which also suffer from small-angle 
approximations breaking down. 
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APPENDIX A: Radial and Transverse Velocity 
Fields 

Given the spherical nature of the problem of redshift 
distortions in an all-sky galaxy survey, it is natural to ex- 
press the velocity field in terms of a spherical vector basis. 
This allows us to set the boundary conditions on the surface 
of the enclosing sphere in an uncomplicated way that avoids 
mixing the 3 cartesian components, and separates out the 
radial components from the transverse components of the 
velocity field for us. In general, we can also decompose the 
velocity field into potential and vortical flows. However, on 
the scales that interest us the generation of vortical flows by 
baroclinic interactions and orbit crossing are negligible. 

The velocity field may be split into radial and transverse 
components by 

r X J , , 

V = Mr , (62) 

r 

where u{r) = r.v is the radial velocity part and J = r x 

V is the transverse term. Given that the velocity field is 
potential, we can define a scalar field, ^'(r), by 

V = V*, 

= r(r.V)* - -(r X L)*, (63) 
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where on the second line we have performed a similar de- 
composition into radial and transverse components of the 
gradient operator, and we have defined L = —ir x V as the 
classical orbital angular momentum operator. 

Expanding $(r) into spherical harmonics we find 

V = ^ Ci„*im„ X 
Imn 

(64) 

where we have introduced the orthonormal vector spheri- 
cal harmonics yJ^^'' and ^[^^ ■ These are defined in terms 
of the scalar spherical harmonics, Yim and obey the same 
eigenfunction conditions; 



L Yim = 1(1+ l)Yim, LzYim = rnYi„ 



(65) 



It is fairly straightforward to show that these can only be 
satisfied by the vectors 



r(L) 



r(M) 



yi(7TT 



r X Y 



(E) 



(66) 



where the factor of 1/^1(1 -|- 1) is required for orthonormal- 
ity. These orthonormality relations are given in Appendix 
C. 



Here we wish to show that the effect of small scale viri- 
alized clusters in redshift space can be expressed by the re- 
lation 



n _ arnm p. 

J^lmn — '^W nn' m' n' ■ 



(67) 



Virialized clusters in redshift space have the effect of smear- 
ing out the galaxy distribution along the line of sight. In 
section 4.3 we have modelled this by randomly assigning 
each galaxy a displaced position, which is then allowed to 
vary like a Gaussian distributed variable. We therefore want 
to relate this picture to the above expression. Perturbing the 
positions of galaxies in the survey by the random displace- 
ment e(r), and averaging over all realisations we find 



Cln 
V 

Cln 
V 

Cln 
V 



APPENDIX C: Orthogonality relations 

When preforming calculations in spherical harmonics, 
it is useful to have reference to the orthogonality, or closure, 
relations harmonic functions obey. In all cases these can be 
expressed succinctly in the notation of Dirac as 



^ ^ (x|lmra)(lmx|x') = 5£)(x — x') 



(69) 



and 



{lmx\yL){yL\l'm'n')ifx = 8^, S'j^^i S^j^^i . 



(70) 



However, in actual application, the symmetry's of the prob- 
lem may require less general notation than that of Dirac. 
Below we have compiled a list of some of the more frequently 
used orthonormal relations. 

The closure relation for spherical Bessel functions is 



ji(kx)ji(kx')k^ dk = — -Sd(x — x'), 



while that for Spherical Harmonics is 



i1m(fi)i^'„/(fi) (iQ = 



(71) 



(72) 



APPENDIX B: Effects of small-scale velocity structures. 'y ^ (Q)YJ* (Q 



The inverse relations for Spherical Harmonics is given by 



47r 



(73) 



where fi is the cosine of the angle subtending Q and Q'. 
Summation over the remaining I modes yields 



J2 Yi^{n)YU^') = 8niA-A') 



(74) 



For systems with azimuthal symmetry these reduce to 



2 

Pi(ji)Pi'(ji) = 21^1 



and 



po(r)S(r)ji(kins')Yi'^(Us')d^r)e, 

po{r}8{r} / p{s - s')ji{kins')Yi"'{flsi)d^ s'd^r, 



+ ll2)Pi(i^)Pi(i^') = 8n(j^- /*')• 



(75) 



(76) 



po{v)8{v) 



I' n' m' 



Swnn'n'iki'n's)YiT (Qs)(i^r, 



Finally, the orthogonality relations for the vector spher- 
ical harmonics of any one type are 



E 



~IV nn' J^l' m' n' , 



(68) 



Yim • Yi*„/rfQ — ^n'^mm'- 



(77) 



where we have expanded the distribution function of the 
e(r)'s as in equation (35). 
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FIGURES 

Figure 1. Coefficients of /5 and Kaiser's Fourier approach 
gives coefficients of 2/3 and 1/5. These are recovered in the hmit 
of A; ^ CO. 

Figure 2. Convolution ^ and redshift distortion matrix V ele- 
ments, showing the mixing of power caused by the selection func- 
tion in the first case and by redshift distortion in the second. The 
figures show elements with / = I': 1=5, n = 3 (solid); / = 10, 
n = 5 (dashed); / = 20, n = 10 (dot-dashed), a) Equal weight- 
ing, convolutions ^^J^^/ (ignoring mask mixing) b) Minimum vari- 
ance weighting (for power) convolutions, showing narrow range 
of wavenumbers sampled, c) Redshift distortion matrix elements 
V^^ I , for /5 = power minimum variance. 

Figure 3. (a) Selection functions for 1.2 Jy IRAS survey, QDOT 
and 2 Jy survey (the last is normalised to unity at small separa- 
tions), (b) weighting functions for 1.2 Jy catalogue: dashed curve is 
minimum variance for the power, assuming (3 = (the 'Feldman, 
Kaiser, Peacock' weighting; dot-dashed and dotted are redshift 
space corrected. 

Figure 4. N-body simulation, analysed in real space, for which 
(3 is effectively zero. The power spectrum is accurately approxi- 
mated by P{k) = {k/ko)~-^-^^ , with A;o = 0.45. The graph shows 
likelihood contours for j3 and ko, picking out the correct values. 
Contour levels are separated by SlnC = 0.5. FIGURE AVAIL- 
ABLE FROM AUTHORS 

Figure 5. N-body, redshift space. Here we assume the correct 
power spectrum, and find the likelihood function for j3. The cor- 
rect value (fi = 1) is within the errors. FIGURE AVAILABLE 
FROM AUTHORS 

Figure 6. Likelihood contours for j3 and the normalisation of the 
real-space number fluctuation spectrum, <7s,IRAS for the IRAS 
1.2 Jy sample. We assume a parametric form for the power spec- 
trum, allowing the amplitude to vary. The data are weighted with 
the minimum variance estimator for power, and modes up to a 
wavenumber k = O.l/iMpc"^ are included. Contour levels are 
separated by SlnC = 0.5. (a) does not include correction for 
small-scale velocity dispersion, (b) assumes an incoherent three- 
dimensional velocity dispersion of 300 kms~^. On these scales, 
the correction is very small. 

Figure 7. Velocity power spectrum, obtained without assuming 
gravitational instability. The shape is, however, assumed, as in 
fig. 4. 
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Selection Function Matrices 




k/(h Mpc ') 
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Redshift Distortion Matrices 
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IRAS 1.2Jy, 2Jy and QDOT Selection Functions 
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